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Abstract: It is not widely appreciated that many subtleties are involved 
in the accurate measurement of intensity-correlated photons; even for 
the original experiments of Hanbury Brown and Twiss (HBT). Using a 
monolithic 4x4 array of single-photon avalanche diodes (SPADs), together 
with an off-chip algorithm for processing streaming data, we investigate the 
difficulties of measuring second-order photon correlations g''^\x' ^t' ^x^t) in 
a wide variety of light fields that exhibit dramatically different correlation 
statistics: a multimode He-Ne laser, an incoherent intensity-modulated 
lamp-light source and a thermal light source. Our off-chip algorithm treats 
multiple photon-arrivals at pixel-array pairs, in any observation interval, 
with photon fluxes limited by detector saturation, in such a way that a 
correctly normalized g^^^ function is guaranteed. The impact of detector 
background correlations between SPAD pixels and afterpulsing effects 
on second-order coherence measurements is discussed. These results 
demonstrate that our monolithic SPAD array enables access to effects that 
are otherwise impossible to measure with stand-alone detectors. 
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Table 1 . Values of first and second order correlation functions for incoherent, coherent and 
thermal light states. Single-mode states are considered, 6 is the angular width of the source, 
A is the wavelength, Tc is the coherence length. Integration effects due to limited detector 
response times and resolution of coincidence counter are not indicated. 

Function g^^^x.z) g^^^x.z) 

Incoherent 1 

Coherent 1 1 



Thermal smc(^7iGj^ exp^-;r^^ 1 ^sinc^ (^nOf^ exp(^-7r^ 

Entangled sinc^^;r0|^ exp^-;r^^, 

[g(2)(o) = l] 



1. Introduction 

The first-order correlation function, g^^\is widely used in optical applications because a simple 
experimental arrangement makes it easily accessible, e.g. Young's double-slit interferometer 
m . In Young's experiment, two points on the wavefront separated by a distance xi2 produce a 
fringe pattern with resultant intensity /i +/2 + 2^Iil2\g^^\xi2,i)\cos{A(pi2) at a screen loca- 
tion where the phase difference is A 912 and the propagation time-difference is T. The fringe visi- 
bility y measures the magnitude of the first-order correlation function \g^^\xi2^T) \ = /^V? ^ 

between field points at the double pinhole. When the light-state distribution, |g(^)(xi2, is 
known then the phase difference is all that is required to define relationship between the field 
points. But in the case of unknown field distributions, measuring this quantity alone is ambigu- 
ous because it cannot distinguish between light states like: entangled photons, incoherent light, 
or coherent and thermal light (See Table [T]). The ability to resolve such potential ambiguities 
is highly significant, for example, for the current debate over the existence of microcavity- 
polariton Bose-Einstein condensation (BEQ.El- lfTOll 

Table [T] shows that to distinguish coherent ifTll and thermal (chaotic) light states |[T2ll 
or entangled photons and incoherent light state, the second-order correlation function 
^(^)(xi2, t) = -^^^^|^^^^^, associated with Hanbury Brown and Twiss (HBT) correlated inten- 
sity fluctuations [tl3J , must also be measured. Here, /i,2(0 the light intensity at a point ±^xi2 
and time t. However, recalling some of the early controversies surrounding photon-correlation 
measurements |[T4l|T5l[T6l[T3, reminds us that second-order correlation measurements are tech- 
nically more difficult to accomplish in the optical range and requires dedicated single-photon 
detectors as well as a coincidence-counting apparatus. An additional difficulty arises due to the 
impossibility of taking single, simultaneous image of the g^^\x^t^x^t) distribution. Instead, a 
pair of detectors and a beam splitter must be used to sample the distribution over a range of 
detector-pair positions. More recent correlation measurements, such as those for BEC cavity 
polaritons, often require the ability to resolve the spatial dependence of g^^^ in a single run. 
Taken together, all these requirements demand an integrated, monolithic, photon detector (as 
flexible as a camera), that is capable of imaging intensity noise-correlations. 

Elsewhere ifTSl . we have presented such a device and applied it in a table-top implementation 
of a stellar HBT interferometer |[T9l , which in their original paper was used to measure the 
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Fig. 1. Micrograph of the 4x4 SPAD array (a), schematics of the SPAD pixel structure (b) 
and electronic readout circuit (c). For a typical operation conditions Vqp = —2\V, which 
by 4 V exceeds the breakdown voltage, Vdd = 3.3y and Vbias = 0^- 

apparent angular diameter of the star Sirius from the current-noise correlations of two detectors 
separated along a variable baseline. As seen from the Table [T] row for a thermal light source, 
the angular width of Sirius star 6 = 0.0063'' assumes that detector separation baseline xn can 
be as large as a 10^ x A, yet enabling to observe second-order correlations. 

However for a table-top implementation, for which the angular width of a typical source is 
6 ~ 1000'' (e.g., a 200 {im spot seen from a distance of 2 cm), the detectors have to be placed 
at a distance xn no more than a few tens of wavelength A in order to observe correlations. It 
was thus not possible to measure g^^"^ in the optical spectral range with stand alone detectors 
located in one plane. A conventional table-top implementation necessitates the introduction of 
a beam splitter together with two detectors at equivalent planes related by reflection, as in the 
original HBT implementation |[T3ll . At the same time for such delicate experiments demand- 
ing a clear-cut answer whether the intensity correlations exist or not, as in the polariton BEC 
experiments, a strong sensitivity to variations in the relative detector position might render the 
measurement procedure impractical and the results biased by an error in the relative detector po- 
sition of the few tens of wavelength scale. Therefore, we proposed L18i the use of single photon 
avalanche detectors (SPADs) integrated into a monolithic array with a simple data-treatment 
algorithm for realizing table-top (laboratory) measurements of an HBT interferometer in the 
original "stellar"-like configuration. 

Here, we report on the peculiarities of our SPAD detector- array for probing local g^^^ correla- 
tions. The structure of this paper is as follows. In Sec.[2l we describe our detector system, which 
comprises a 4x4 array of Si SPADs implemented in CMOS technology (Fig.[T]) and a simple 
algorithm to treat the multiphoton time-of-arrival distributions from different SPAD pairs, im- 
plemented off-chip. The device also incorporates on-chip high-bandwidth I/O circuitry which 
facilitates the external data processing. We discuss the data treatment procedure to correct for 
spurious correlations in measurements (Sec. [3]) and we demonstrate the operation of our device 
for cases of very different light statistics (SecO: measuring the beat note of a multimode co- 
herent source, measuring the depth of intensity oscillations in a superposition of incoherent and 
thermal source, and measuring the instrument response function of a table top interferometer 
as a function of the source angular- width. 

2. Structure of g(^) imager 

Our g(^) imager is based on a single photon avalanche diode (SPAD) array ll20ll . whose pho- 
tomicrograph is shown in Fig. [T](a). The array pitches are 30 and 43/im in the horizontal and 
vertical directions, respectively. The chip is fabricated in a 0.35 /im CMOS technology. Each 
pixel of 4x4 SPAD array is independently accessible and has a structure as indicated in Fig. 
[T](b). The SPADs p+-n- well junction is reverse biased above breakdown by a voltage known 



as excess bias voltage. In this mode of operation, called Geiger mode, single photons may be 
detected and counted using auxiliary electronics. The anode is enclosed in the p+ -diffusion ob- 
tained from the standard implant used to define a PMOS channel. To prevent premature edge 
breakdown in the periphery of the multiplication region, a 3.5 /im guard ring of lightly doped 
p-type implant is designed to surround the anode. 

The p+ anode of each SPAD pixel is kept at a negative voltage Vop of 21 V, which by 
4 V exceeds the breakdown voltage (Fig. \T\ (c)). When an avalanche is triggered, e.g. as a 
result of photoabsorption in the depleted n-well region or thermally induced carrier release 
from a trap defect, the avalanche is subsequently quenched via a resistive current path in the 
PMOS transistor structure with the transistor gate kept at a ground potential ( Vrias = 0)- This 
ballast resistivity brings the bias voltage across the p+-n junction below the breakdown voltage 
and is used to read out the photodetection events. After avalanche quenching followed by a 
recharge of the depletion region, the detection cycle is completed and can be initiated again by 
arriving photons or thermally excited carriers. Each SPAD pixel comprises a high-bandwidth 
inverter so as the Geiger pulses are converted into digital signals and transmitted for off-chip 
data processing. This approach considerably improves the signal-to-noise ratio, however as 
discussed in section [51 the data should be corrected for spurious correlations due to a possible 
crosstalk between SPAD pixel outputs. 

During the time required to complete a detection cycle, the detector is insensitive to arriving 
photons. In our structures, the dead time is 12-15 ns, so as a theoretical limit on the count rate 
set by detector saturation (pile-up) effects is about 30 MHz. 

Tunneling of minority carriers through the shallow p-doped guard ring and/or trap-assisted 
carrier release may trigger an avalanche and thus create spurious Geiger pulses. The lowest 
detectable photon flux is set by the dark count rate (DCR) of SPADs, which at room temperature 
conditions, is in the 5-10 Hz range. Such low DCR is achieved by implementing n- wells of 
small diameter d=3.5 /im. The lowest detectable photon flux density in our experiments is thus 
10~^ photons/s/cm^ and the dynamic response range is 64 dB. 

At 4V excess bias voltage, the measured photon detection probability (PDP) of such 0.35/im- 
CMOS SPAD pixels is 40 % at a wavelength 450 nm and exceeds 25 % in the wavelength range 
at 400 - 550 nm. PDP is proportional to excess bias. The spectral sensitivity curve and pixel 
statistics on dark count rate can be found in [20j. For our experiments, the chip was mounted on 
a standard ceramic case and was wire-bounded with 20 /im-thick aluminum wires of 6-7 mm 
length spaced by 100 /im. The bonding wires can be seen on the peripheries of the chip in Fig. 
Fig.[T](a). For what follows it is important to know that the bonding wires of individual detector 
data lines are ordered with the SPAD pixel number / indicated on the microphotographic image. 

Since all 16 detectors in the array have separate parallel outputs, (^2^) = 120 HBT interfer- 
ometer measurements are possible simultaneously between various detector pairs at temporal 
resolution set by the SPAD timing jitter characteristics ( ^Ops). Measurements of pairwise sec- 
ond order correlation are based on implementing expression II2TII : 

(2)/ ^ Tr{pa^ajajai) 

TY[paiai)TY[pajaj) 

with a simple algorithm. Here, integers ij {i ^ j) enumerate detector pixels Di and Dj in 
Fig. [U dij is the photon annihilation operator at two detectors, p is the density operator for 
the field and detectors, and trace is taken over the detector and field states [22] . At this stage, 
the pairwise intensity noise correlations g^^\xij^T) are computed from multiphoton arrivals 
at different detector pairs of the array, using a programmable four-channel 6-GHz bandwidth 



digital oscilloscope (Wavemaster 8600A, LeCroy) and implementing a simple expression 



Lm=0 Lj 




(2) 



where Xf and Xj are discrete random variables whose values (no event) or 1 (photon detec- 
tion) correspond to the binary data stream emanating from any pair of detectors Di and Dj, 
respectively. The spatial lag Xij is set by the separation of the detector pair within the SPAD 
array. Time-lag increments z = IT are set by multiples of temporal resolution T (equivalent to 
the bin width in conventional methods), where NT is the width of the measurement window 
and M is the overall number of measurements series. Note that both addition and the bitwise 
AND operator (A in the numerator of Eq.[2l) could be implemented easily using on-chip digital 
electronic circuits. 

As can be seen from Eq. [2] , the concept of our g^^^ imager ifTSl is drastically different to 
conventional approaches based on measuring start- stop times and computing a timing his- 
togram of delayed photon arrivals at two detectors, followed by renormalization of the his- 
togram according to a statistical hypothesis for the field (e.g., a hypothesis like g^^\T = 0) = 1 
or g(^) (t = ±NT /2) = 1). Our g^^^ -imager algorithm (Eq Q) does not suffer from missing de- 
tection events in case of subsequent photon arrivals at the same detector and operates correctly 
with multiphoton arrivals. Here, we will show that differently from the conventional approach, 
it operates correctly with modulated intensity signals at any count rates within the dynamical 
response range of SPADs and any width NT of the temporal window of interest. 

An inherent limitation of our current prototype detector- array is the off-chip data treatment 
implemented with a general purpose digital oscilloscope that enables user-defined Matlab data 
processing subroutines. As a result of the relatively slow Matlab interpreter, a data trace of 5 
lis length (N=5000 points at T=lns resolution) requires 75 ms of computations to calculate 
the second-order correlation function in 200 points (in the interval —100^^ < T < 100^^). To 
achieve a low standard variation of the processed f data, about M=10^ such data traces needs 
to be acquired and treated. The undesirable latency of this numerical processing step could be 
reducesd dramatically by integrating our robust algorithm into the chip. 

3. Multiphoton correlation measurements 

The performance of our g^^^ imaging device (Fig.[T]) was tested by measuring the distributions 
g^^\xij^l) when the 4x4 SPAD array detector is uniformly illuminated with an incandescent 
light bulb (incoherent broadband light source). Imperfections in g^^^ imaging with the SPAD 
array chip might be caused by several reasons including the concept of monolithic implemen- 
tation itself. Fig. [2] shows two time-resolved correlation patterns g\j ^ (t) (left panels) and corre- 

lation maxima images g\j ' (0) (right panels) acquired with respect to detectors D8 (top panels 
, J = 8) and D9 (bottom panels, j = 9) respectively. In the figure, the pairwise correlations 
g |g ^ (t) {g^^ (t)) are ordered according to the detector index / in Figdl For an ideal g^^^ -imager, 
one would expect to measure no correlations {g^^\xij^T) = 1) and a uniform correlation max- 
ima map for all pairwise detector combinations. However, Fig. [2l reveals a deviation from the 
expected distribution. 

Measured second-order auto-correlation function of a detector f ^(t) shows a large coin- 
cidence peak g^^\0) 1 at zero time lag [green scale curves in (a)]. In fact, measurements 
with one and the same detector at T = do not represent statistical correlations between field 
states projected on two independent detectors in Eq. [T] (211 l22l- Instead, the statistical algo- 
rithm in Eq.[2]is applied to one and the same field state projection. The excess correlations at 



D12&D8 




(b) 



-20 -10 10 20 -20 -10 10 20 -20 -10 10 20 -20 -10 

(C) 



Time [ns] 
10 20^ 



DO D4 



D1 D5 



D2 D6 D10 D14 





DO D4 




D12 


D1 


D5 




D13 


D2 


D6 




D14 


D3 


D7 


D11 


D15 



Fig. 2. Second-order correlation functions g^^\xij,T) of SPAD array pixels (a) and the 
maps of correlation maxima g^^\xij,0) [(b) and (c)] measured for incoherent light using 
detector D8 [(a) and (b)] or D9 (c) as a reference. The time resolution is 1 ns. The panels 
corresponding to individual SPAD pixels are arranged in the same order in which they 
appear in the imager pattern and corresponds to SPAD pixel position in Fig[T] Green scale 
curves in (a) show autocorrelation curves g^^\0,z) in linear (olive curve, left axis) and 
logarithmic (green curve, right axis) scales 



zero time-lag can be understood by noting that the bitwise AND-ing of the binary data- streams 
in Eq.0 yields {X^\n) AX^'^\n)) = {x\^\n)) = jlT for the average detection rate jl and 
time window NT. The result is then g\f\o) = ^^^^^j =l/jlTAn FigO the correlation curves 
are acquired at resolution time T=l ns and the average rate of detections 6 MHz, yielding an 
estimate gfj{0) = 166, in perfect agreement with the acquired data in Fig 121 (olive curve, left 
axis). 

(2) 

In the vicinity of the self-correlation peak, the dead-time bands of g\- \t) = are due to a 
recharge of the detector depletion region after each Geiger pulse. Measurements in our SPAD 
array indicate the dead time T/) of 12-15 ns (Fig|2l(a), green curve, right axis). At larger time 
lags |t| > Tz), the self-correlation function g\f\T) again shows an excess of delayed coinci- 
dences, exponentially decaying to g*^^) = 1. [Note a logarithmic scale of the right axis in (a).] 
This time, i.e., the excess of delayed coincidences, is due to the thermal release of minority 
carriers from traps in the depletion region, that were captured during the preceding avalanche 
current pulse. A casuality relationship between afterpulses and their generating photodetection 
events explains the measured increase of delayed coincidences for the incoherent light source 
EH EH: 

where PA{t)dt is the probability to detect an afterpulse in the time interval {t^t -\-dt) and £ = 
J^^PA{t)dt is the overall probability of afterpulsing. We stress that the impact of afterpulsing 
effects on auto-correlation function ^ varies with the incident photon flux and is particularly 
pronounced at small count rates fl. Measuring auto-correlation curve ^ for a larger time lag 
(not shown in the figures) shows the characteristic decay time of afterpulsing probability 40 ns 
with an integral effect e = 7% at 4 V excess above threshold. The afterpulsing probability e 
is determined by the number of trap defects in the active (n-well) region material, the number 
of carriers involved in an avalanche and the recharge time. Its impact on ^ can be drastically 
reduced by increasing the deadtime of the detector, e.g. by using a smaller recharge current for 
the depletion region. We have also noticed a slight decrease of afterpulsing probability with 
excess of the biasing voltage above the avalanche breakdown. In any case, the auto-correlation 
measurements, as in ([3]), do not allow the second-order correlation statistics of a source to be 
accessed but can be used to measure the count rate, dead-time, and afterpulsing probability of 
SPAD detectors. 

The results of pairwise measurements in Fig|2l deviate slightly from the expected cross- 

(2) 

correlation function gfj (t) = 1 as well. For / ^ j, the axis scale is limited to the range 

(2) 

0.8 < gij^ < 1.3. In each of the correlation maps in Figs, (b) and (d), one can distinguish dis- 
crepancies of two kinds: (i) There are a couple of SPAD pixel-pairs exhibiting a bunching of 
Geiger pulses at zero time-lag (red curves), while (ii) the rest of detector pairs exhibit anti- 
correlations at T = (black curves). 

Thus in Fig|2l(a), which was acquired using SPAD pixel D8 as a reference, detectors D7 and 
D9 exhibit correlation maxima g^^'' (0) ~ 1 .3 measured at temporal resolution r=lns. At higher 
temporal resolution and shorter integration time (r=100ps), its height increases while the nat- 
ural FWHM measures 320 ps in width (not shown in the figures). Measuring these spurious 
correlations at various excess voltages above the avalanche threshold and incident light inten- 
sities, we have not seen any impact. Such behavior attests for the electrical cross-talk origin of 
spurious correlations that occur at the SPAD output data lines and have no impact on the detec- 
tion process. In particular, the detectors D7 and D9 have data line bonding pads on the CMOS 
chip adjacent to the wire of SPAD pixel D8. (The bonding wires of SPAD data outputs seen at 
the chip edges in Fig [T] are ordered with the SPAD pixel index /.) These wire lines are subjected 



to high currents when the logical state changes between (0 V) and 1 (3.3 V). The measured 
excess of correlations can thus be attributed to electrical cross-talk between the SPADS pixels 
with neighboring bonding wires. This is confirmed in measurements of pairwise correlations 
with a reference detector D9 (Fig.[2](c)) showing the excess of correlations for detector pairs 
D8-D9 and D9-D10. The patterns in Fig 12] indicate that there is no direct optical crosstalk ll25l 
between the adjacent SPADs within the array, which would appear in Figs|2](b) and (d) as a 
cross pattern centered at the reference pixel. 

All other SPAD pixel pairs in Fig |2l with larger index difference (| / — j| > 2) have non-nearest 
neighbor bonding wires and show small anti-correlations (0) ~ 0.85. The anti-bunching dip 
of Geiger pulses is of 2 ns width followed by strongly damped oscillations of 4 ns period 
(they can be seen in Figl3] (b)). This feature corresponds to reflections in the 90 cm length 
RF cables transmitting photodetection events data to the processing timing electronics (oscillo- 
scope Wavemaster 8600A, LeCroy). In particular, we attribute this to the impedance mismatch 
between CMOS output inverters of the SPAD array chip and the 50 Ohm RF cable. Measure- 
ments of spurious anti-bunching dips at various illumination power and excess above threshold 
conditions have not revealed any variations in the dip height or width, attesting the non-detector 
origin of measured g JJ^'^^(t) bias. 

Spurious correlations and afterpulsing effects have to be taken into account in the photon 
coincidence measurements utilizing our SPAD arrays. Starting from the results of Ref. 1231 and 
spurious background correlations g\f'^^{'^) between Xi and Xj data lines, one can conclude that 
for incoherent and coherent light sources as well as for thermal, entangled or single photon 
sources with coherence time smaller than the detector dead time Td, the measured correlation 
function (O is 

= 1 + tttV ^ ^'^'"'^^ " ^"^^ 

(2) 

where g^j^ is the correlation function for photons that would be measured by ideal detectors. 
The second term takes a reduction of the photon bunching peak (anti-bunching dip) due to the 
spurious afterpulses of one detector uncorrected with photon detection events and afterpulses 
of another detector. Note that spurious uncorrected dark counts have the same effect (see also 
Eql5]), however for illumination intensity levels used in the experiments discussed here, their 
impact is negligible. The third term is the correction due to background correlation caused by 
data transmission lines from the SPAD chip to the timing electronics. 

We shall stress that in agreement with Eq.(|4|, there were no variations observed in the impact 
of afterpulsing effects on the cross-correlation function with increasing (decreasing) incident 
light power. This feature allows the detector's pairwise correlation measurements (O to be 
used as an efficient tool for accessing the field statistics g^^\xij^T). Unlike conventional de- 
tection methods based on start- stop timing histograms of delayed single photon arrivals, which 
are not capable to treat multiple photon arrivals at one detector during the start- stop inter- 
val measurements, our approach implements properly normalized multiphoton distribution. As 
such it is robust against missing detection events, the impact of Poisson-like distribution decay 
oc exp(— /ir) (see Ref ifTSll for details) and intensity modulation. As opposed to standard meth- 
ods, Eq.Opermits arbitrary count rates and temporal window of interest and does not require a 
statistical hypothesis to normalize g''^\ 

4. Experimental results and discussion 

The imager was tested by measuring statistical properties of several light states including mul- 
timode coherent state, near field of an extended intensity-modulated thermal light source and 
the far field of stationary thermal light source. 



4.1. Multimode coherent state 



As a model system for a multimode coherent state we use the emission of a He-Ne laser 
operating on a fundamental Gaussian mode of the cavity at 633 nm wavelength transition. 
The cavity length of L=2l cm assumes a large separation of the longitudinal modes c/2L=710 
MHz, such that a small optical gain and narrow linewidth allows only two or three longitudinal 
modes to reach the lasing threshold (in function of detuning from the gain line center). No 
special provision for stabilization of the cavity length or lasing modes has been made. The 
relative phases of modes vary rapidly in time such that no intensity modulation at the intermodal 
frequency can be observed in the emitted output radiation. 

The output beam of the laser impinging the g^^^ imaging detector was attenuated to reduce 
the count rate at detectors down to 2 MHz. The temporal resolution of the timing electronics 
was set to r=100 ps so as to observe dynamics of lasing modes making a cavity roundtrip in 
1.4 ns. 

Figl3](a) shows the correlation function gf\{"c) [Eq.(|2|)] measured by two detectors in the 
middle of array (detectors D5 and D9) at 30 jim baseline. The background component due 
to a spurious crosstalk at the same detector pair ^^^^'^^(t), which was measured with an in- 
coherent broadband light source, is shown in (b). The corresponding correlation function of 
the multimode laser beam g^^\xs^9^ t) calculated from Eq.© after corrections for background 
correlations and afterpulsing effects in the detectors is plotted in (c). It can be seen that this 
procedure is an efficient tool for removing spurious (anti-) correlations seen as a dip in both the 
measured correlations gff'^'c) and reference background f [^'^^(t). 

According to the Table [T] for each of the modes of the laser, one measures g^^\T) = 1. 
However, as evidenced by the measured second order correlation function in Fig|3](c), when 
these modes are simultaneously excited in the lasing spectrum, one can see the effective beat 
signal of intensity and phase noise correlations in the modes |[3Qll . The numerical fit with a 
cosine function (green curve) reveals the noise beat period of 1.4 ns, that is the time of a 
roundtrip in the laser cavity. The amplitude of harmonic oscillations 0.2 indicates that the noise 
of different modes is just partially correlated. We stress that no such periodic signal is observed 
in the intensity of the output laser beam. 

4. 2. Intensity -modulated extended thermal light source 

As a model system for an extended quasi-monochromatic thermal (chaotic) light source, the 
green line of mercury (546 nm) from the emission spectrum of a low-pressure Hg-Ar discharge 
lamp (CAL2000, Ocean Optics) is used. The lamp is driven with an original power supply 
incorporating a 30-KHz AC source producing 1.8 KV pulses to initiate the discharge and 120V 
in a steady state operation conditions. The green line of mercury is filtered with a 10 nm band 
pass filter (FL543.5- 10, Thorlabs). 

Such a quasimonochromatic chaotic light source exhibits the first order correlations \g^^^\ > 
in the Young's interference experiment (Fig|4l(a)). It is characterized by Gaussian distribution 
of photons in the Glauber /^-representation (Plank distribution in ^-representation) so as one 
should expect measuring second order correlations g^^^ (x, t) = 1 + ^ \g^^^ (x, t) p with q being 
the number of photon modes (of equal intensities) reaching the detectors ifTTIl . Here we report 
on measurements in the near field of the source such that the number of modes is large. (The 
angular width of the source ~1 in Table [T]). Therefore, for a detector pair with the baseline 
of the order of A and more, that is practically in all cases, the second order correlation function 
should be close to one, indicating no second order coherence. 

To measure the second order correlations in the near field of such chaotic light source, a 3 
mm size domain of the lamp discharge region was imaged onto the SPAD array at magnification 
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Fig. 3. Multimode coherent state emitted by a He-Ne laser at 632.8 nm wavelength, (a): Ac- 
quired second order coherence function g\g{"c) [Eq.(|2])]. (b): Spurious correlations back- 
ground gf'g^^i'c) measured with the help of an incandescent light bulb, (c): Second order 

correlation function of the field g*^^)(x5 9,T) corrected for spurious correlations [Eq.(|4])]. 
The laser cavity length is L=2l cm with the cavity roundtrip time of 1.4 ns (the period of 
oscillation seen in (c)). The detectors D5 and D9 separated by 9=30 jUm baseline are 
used at a temporal resolution of timing electronics r=100ps. 





Time, [is 

Fig. 4. Extended quasi-monochromatic thermal light source (546nm line of mercury), (a): 
Young's interference fringes indicate phase correlation > 0. (b) and (c): Measured 

(2) 

second order correlation function g in the near field of the source obtained after correc- 
tions © at temporal resolution lOOps (b) and 100ns (c). The oscillations are due to the AC 
power supply of the Hg-Ar discharge lamp. 



0.44. The image created by a 150 mm focal distance lens was overilluminating the entire SPAD 
array such that we were able to test correlations at two pairs of detectors, one in the middle 
of the array (SPAD pixels D5 and D9) and one pair with detectors at the opposite array cor- 
ners (pixels DO and D15). These detector pairs with baselines of 30 /im and 158 /im measured 
correlations t) at the points of the lamp discharge that are separated by 68 and 360 /im 

distance, respectively. We were expecting to see no correlations (g^^ = 1), since the correlation 
measurements were done between different emission points of the extended source. However, 
in both cases, we have observed excess of coincidence counts. FigHdetails the correlation func- 
tion ^5 9(t) measured at the detector pair in the array center. The correlations measured with 
SPAD pixels at the array corners exhibit the same correlation curve with excess of coincidences 

^0^15 (^) ^ ^^^^ shown in the figure). 

The difference between the usual statistical data analysis and our method, can be under- 
stood as follows. The distribution g^^\T) is assumed to be an even function about the origin 
with boundary condition g^^\T ±oo) 1. However, when applied to a given measurement 
window of interest that is at T = ±NT /2, this assumption can be invalid. Consider the high 
resolution data in Fig. |4](b) sampled at r=100ps. Those data, which may be regarded as a 
uniformly distributed histogram of two-photon arrival times, can be renormalized to meet the 
above boundary condition, viz., g^^^t) = 1 at the edges of the sample window. But because of 
the renormalization, will not have the required value of 2.7. 

By limiting the temporal resolution to 7=100 ns and increasing the width of the measurement 
window to NT=50 /is, the origin of the excessive correlations can be made visible [Fig.|4](c)]. 
The measured second order correlation function just reveals the fact of long-period intensity 
oscillations of the lamp at the double frequency of the AC power supply, at about 60 KHz. 
Because of these oscillations and the large angular width of the source, the photon bunching at 
T = due to Planck's statics 1 13 ] is not visible. The change in the width of the measurement 
window and the temporal resolution (in Figs, (b) and (c) ) has no impact on the measured 
excess value g^^\0) = 2.7, attesting for correct normalization of the correlation function when 
multi-photon arrivals are treated with the algorithm (O-©. Interestingly, after correction for 
the afterpulsing effects at the detectors [£ in Eq.(|4|)], the minima of g*^^^(T) function are close 
to zero, in agreement with fluorescence lifetime considerations. 

Modulation depth of the correlation function can be tailored in a mixed state. Fig. [5] details 
second-order correlation measurements for a superposition of incoherent light (as in Fig 121) 
and a quasimonochromatic chaotic light beam (as in Fig I?]) at various probabilities to detect 
photons in chaotic light state. The incandescent lamp is used as a source of incoherent emission, 
which is superimposed on mercury green line fluorescence from the intensity-modulated Hg- 
Ar lamp. The probability P{Hg) that a detected photon has been emitted by mercury atoms, is 
given by the intensity ratio lug/l estimated from the count rates at the detectors. For example, 
P{Hg) = 0.38 when lug = 6.8 kHz, / = 11 MHz. 

One can easily show that for a superposition of mutually uncorrected photonic states, char- 
acterized by partial intensities I^^^ and correlations ^(^^'^(x/j, t), the second order coherence 

function g^^\xij,i:) = ^j^^^j^^^ for the overaU photon fluxes lij = Ea^y ^ impinging the de- 
tectors assumes the expansion : 

This expression can be simplified in the case of equal count rates of the detectors, yielding 
the partial weights (/(^))^/(/)^ for state contributions of ^(^^'^(t) — 1 to the overall excess 
(or lack) of coincidences of the correlation function. For example, for the overall rate (/) = 



(a) (b) 

3l ^1 




• ^ • ^ • ^ • ^ 0-^ , ^ , , 

-2 -1 1 2 0.5 1 

Time x, |is Squared relative intensity, /^g //^ 

Fig. 5. Intensity correlations measured for the superposition of incoherent hght (from 
incandescent lamp) and intensity-modulated quasimonochromatic thermal light (green 
546nm line of mercury from Hg-Ar lamp CAL2000) at various relative intensities. Tempo- 
ral resolution is r=100ns. (a): Second order correlation function g^^^T). (b): Correlation 
function maximum g*^^) (0) plotted as squared relative intensity of thermal light incident on 
detector pair. 



(1 +e)(/(^)) of photodetection events caused by a process of intensity and statistics 

g\f'^{^) accompanied by spurious (incoherent) afterpulses occurring with probability e, Eq.© 
predicts a reduction of the bunching peak (or antibunching dip). The second order coherence 
function of the entire process is then 1 + {g\f'^{^) ~ + ^)^^ in agreement with EqlH 
Along the same lines, the effect of spurious dark counts of detectors can be taken into account 
at low count rates (lij). 

For the experimental results reported in FigE the incoherent light constituent does not con- 
tribute to the excess of correlation (the correlation function g(^)^^^'^^^(^z) = 1). The average inten- 
sities at the detectors (SPAD pixels D5 and D9) are equal so as the contribution from the chaotic 
light source of Fig|4]is reduced as 1 + {gf'g^^{T) — l){lHg)^/ {I)^ [Fig|5](a)]. The maximum 
values of the second order correlations ^ (0) — 1 linearly grows with the squared probability 
of detecting a thermal light photons {Ing)^/ {I)^ [Fig IS] (b)]. 

4.3. Quasi-monochromatic thermal light source 

So far we have not yet observed the HBT photon bunching in the emission from a quasi- 
monochromatic chaotic light source. To make it visible, we replace the power source of the 
lamp and conduct measurements in the far field. In particular, we drive the Hg-Ar lamp bulb (the 
same as used in previous experiments) from a DC voltage source. This lamp bulb (CAL-2000- 
B, Ocean Optics) having a cold cathode and U-folded discharge and designed for operation 
in the AC regime, also operated well with a DC power supply (160V @ 15 mA). To start the 
discharge, we used the original AC power supply of the lamp, which was connected in parallel 
with a DC source via a filter (2H inductance) such that the AC supply was gradually turned off, 
whereas the DC source was gradually turned on. 

As in the experiment of FiglH the light emitted by the lamp was transmitted through a lOnm- 
width bandpass filter, which keeps only the emission at the green line of mercury (546 nm). 
The inherent difficulty is related to the fact that the number q of uncorrected photon modes 
contributing to the intensity interference at SPAD array detectors g^^^ (•^, t) = 1 + ^ (x, t) p 
has to be kept small. (The requirement q \ can be also seen from Eq.© for a superposition 
of q modes of equal intensities.) For a light beam of transversal size w seen from detectors at a 
distance L in the ray cone of angle = xn/L, the number of photon modes is ^ = / ~ 
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Fig. 6. (a) Experimental setup of the table- top stellar HBT interferometer, (b) Correlations 
measured at various detector separation Xij and distance L to the fiber end. The Fresnel pa- 
rameter FN = is indicated in the panels, (c) Measured (points) and modeled (curves) 
second-order correlations in function of detector separation for the model Eq|6l(dashed blue 
curve) and Eq. [8] (solid red curve). The inset shows the corresponding near field distribu- 
tions at the fiber end (indicated with same type ). (d) shows the corresponding modelled 
g(2)„,ax^j,^ patterns (top line of panels) and far field patterns (bottom line of panels), (e) 
Imaged second-order correlation maxima along the row of the array. Its position in the 
g*^^) plane is indicated in (d) with green solid lines. It is assumed that g^^'^ = 2 along the 
diagonal. The green line of Hg (546nm) is used. Temporal resolution T = 1 ns. 



^{kw(^f = Tei^^)^- For a SPAD array of a few 100 x 100 jim^ size, which defines a 
characteristic basehne xi2 of detectors, the requirement q ^ 1 imposes the angular width of 
such extended source to be = w/L ~ 10~^ or less (see also Table [T]). 

To compromise the brightness of the source and its angular width at the detectors, the light 
from the lamp bulb was injected into a multimode fiber of 1 m length and core diameter w ^ 200 
/im [FigHa) ]. The other end of the fiber was used to illuminate the SPAD array in the far field 
zone of the fiber end, at a variable distance L (1 - 3 cm) from the fiber such that the whole 
4x4 SPAD array is over illuminated. Such extended thermal light source is of the small angular 
width w/L= 10~^ rad and exhibits first-order correlations g^^\ when the SPAD array is replaced 
by Young double-pinhole interferometer as in Fig. [4] (a), which was taken in the same setup 
configuration but with the Hg-Ar lamp driven from the AC power supply. 

In intensity interference measurements with a quasi monochromatic thermal light source, the 
coherence width and Doppler-broadened spectral line have the same impact as in the case of 
the Young double slit interferometer. The classical expression for second-order spatio-temporal 
correlations for a non-polarized single-mode chaotic light source are determined by the first- 
order correlations |[TT] |2T| |26l with coherence time Tc=2^/2n\n2 / Ao) due to the inhomoge- 
neous broadening Aco (FWHM) of the emission line 



where the second term in the right hand side takes into account the decorrelation effects due to 
unpolarized light (the coefficient 1/2), the zero-delay degree of spatial coherence and the Gaus- 
sian profile of the delayed first-order correlation function for an inhomogeneously-broadened 
line. 

FigEl (b) shows several traces g^^\T) measured at different distances to the fiber end at 
several detector pairs (the Fresnel like parameter FN = ^ is indicated for each trace). At 
small FN < 0.5 (large distance to detectors and small baseline), the detector counts show a 
correlation peak. At FN = 1, it disappears and then at FN ~ 1 .2, it reappears but with a smaller 
amplitude of correlation excess. At first sight, such behavior corresponds well to the expected 
oscillations of the sine function in Eq.®. However measuring correlations at various distances 
from the fiber end with several detector pairs and correcting the data for afterpulsing effect 
(FiglSl (c), points), we find large discrepancies with expected dependence (Fig|6l (c), dashed 

(2) 

blue curve). The difference is outside the instrumental error and the fast growth of gmax at 
FN ~ 1 .2 cannot be compensated for with a numerical fit to the model in Eq.®. 

To interpret the experimental results in Fig|6l one shall admit that g{2), being a function of 
two variables, can be tailored not only in the time domain, as exemplified in previous sections, 
but it can also be tailored via spatial modulation of the near field distribution. This fact is often 
omitted in experiments on photon bunching, leading to wrong data interpretation. In the same 
way as the time domain variations of intensity (intensity noise) renders invalid the exponential 
term of the famous expression the sine term in that expression is obsolete if the source has 
a nonuniform intensity distribution in the near field. 

Using Glauber' expression g^^^ = 1 + g[^^ G^2i / g[\^ G^^ O with g|]^ = {EiE]) for 
an extended source with near field intensity distribution I{x,y) and detector pair with 
baseline Xij centered at x = 0, after integration over the wavelet contributions £"1^2 
j^\Eo{x,y)\^Qxp{—ikxi2x/L)dxdy from mutually incoherent elements dxdy of the extended 
source, one finds 




(6) 




(7) 



which is nothing else than the Van Cittert-Zernike theorem applied for visibility of the sec- 
ond order correlations. Here, as in ©, the coefficient 1/2 takes the decorrelation effects due 
to unpolarized light into account. We shall stress that only the 1-D Fourier transform of the 
near field intensity in the direction of the detector baseline is involved. For a uniform intensity 
distribution, integration over the source area from — w/2tow/2 yields expression (O. 

Field distribution in multimode fiber highly depends on excitation conditions and uniformity 
of the fiber core material 1 27 , 28 . 29 1 . Due to cylindrical symmetry of the fiber, the intensity dis- 
tribution has an intensity bump or a dip at the core center. To simplify analysis, we consider only 
the main radial harmonic of frequency Q., assuming the approximation /(r) oc 1 + Ycos{Q.r): 



g(2)(x,,-,T) = l + i 



smc(^nFN^ + |sinc(^7rFA^+ ^Q.w^ + |sinc(^;rFA^- ^^w^ 



exp 



(-4) 



1 + 7sinc(^^lw) 

(8) 

where yis the intensity modulation depth and we have used Fresnel parameter FN = wxij/XL 
to shorten the expression. The solid red curve in Figl6l(c) shows the result of the numerical fit 
with the model dU), reporting 7= —0.3 and spatial frequency ^1 = 2.8 ;r/w. The corresponding 
near field (shown in the inset) has a small intensity dip at the center. Note that the difference to 
uniform intensity distribution is remarkable only in the correlation maxima distribution gmax{^) 
(Figl6](d), top line panels) while it is not visible at all in the far field distribution impinging the 
detectors (Figl6l(d), bottom line panels). 

Being limited by the number of acquisition channels, we were able to record simultaneous 
correlations between four independent detectors form one row of the SPAD array. In Fig|6] 
(e), g^^\xij^Q) measured along the array row is plotted as a pairwise correlation map g^^\ij) 
for measured (left panel) and calculated with the models ([8]) and © correlation maxima (two 
right panels). In this image map, the spatial oscillations of the coherence factor are clearly 
visible. The measured excess of correlations for the corner pixels on the secondary diagonal 
corresponds better to the model ([5]), confirming previous conclusion on impact of the spatial 
modulation of the light source. 

5. Conclusion 

We have presented a g(2)-imager built with conventional CMOS technology, which is capa- 
ble of measuring second-order spatio-temporal correlated photons. We have discussed several 
important aspects related to temporal and spatial modulation of the source intensity when con- 
ducting HBT correlation measurements with such imager. Such detectors will find various ap- 
plications. This approach allows the functions g^^^^ of other even orders to be implemented as 
well. 

The application in Sec. 14.31 used a green line emission of mercury in a Hg-Ar discharge 
lamp as a reliable source of HBT photons pairs. The lamp was emitting around hundred /iW 
of thermal light at the green-line transition, only small fraction of which can be coupled into 
a multimode fiber. Yet the photon flux ji emanating from the fiber was quite high, about 10^^ 
photons/sec. If not decorrelations due to the angular width of such source, one would need 
detectors and coincidence electronics with integration time T ~10ps to detect coincidences at 
the highest rate Rc ^ Iji^T of about 10^^ photon pairs/sec in the near field of the fiber. 

In practice, because of the low collection efficiency of detectors in the far field of the fiber 
end, the count and coincidence rates are very low, about /i=100 kHz and Rc=\Q Hz respec- 
tively, for integration time T=\ ns. Using detectors of large diameter d will not overcome the 
problem, yielding a reduction of measured excess of correlations {g^^\"c) — 1) by a factor of 
~ w^J^/A^L^ [see Sec l4.3l 1, like long integration time T results in a decoherence factor of 



Although the alternative of building a quasi-thermal source by using a He-Ne laser beam im- 
pinging a rotating sintered disk offers higher intensity, larger coherence time, and hence much 
faster acquisition, we scrupulously avoided this type of source so as not to introduce the un- 
desirable possibility of spurious correlations such that to observe a clear impact of detector 
imperfections, non-stationarity and non-uniformity of the field on measured second-order cor- 
relations. 

It is interesting to stress the impossibility of detecting a normalized correlation function 
g(^) (t) for entangled biphotons (Table[T]) |[T2l . A state of the art entangled photon source based 
on spontaneous parametric down conversion in periodically-poled nonlinear crystal shows 
conversion efficiency of 10~^ and is capable of producing bi-photons at rates of 10^ photon 
pairs/sec per miliwatt of pump power. The principal difference to a thermal light source illumi- 
nating detectors is that all generated bi-photon pairs produce coincidence detections, provided 
perfect collection efficiency and no decoherence due to the angular width of the source. Tech- 
nical difficulties arise from low overall power and a very short coherence time of bi-photons 
(typically, ~0.1 ps) so as to measure the peak value g^^^O) = 1 for entangled photons (see 
Table [TJ one will need detectors with the response time T ~0.1 ps, which do not exist yet. On 
practice, one just measures non-normalized second order correlations G*^^) (t) with the time lag 
T defined by the optical path difference to two detectors. 

Future work will include the development of larger arrays of SPADs, the integration of on- 
chip data processing based on equation (1), and the extension to other g^^^ -imaging applications. 
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